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ML-PDA  and  ML-PMHT:  Comparing  Multistatic 
Sonar  Trackers  for  VLO  Targets  Using  a  New 
Multitarget  Implementation 

Steven  Schoenecker,  Peter  Willett,  Fellow,  IEEE,  and  Yaakov  Bar-Shalom,  Fellow,  IEEE 


Abstract — The  maximum-likelihood  probabilistic  data  associa¬ 
tion  (ML-PDA)  tracker  and  the  maximum-likelihood  probabilistic 
multihypothesis  (ML-PMHT)  tracker  are  tested  in  their  capacity 
as  algorithms  for  very  low  observable  (VLO)  targets  (meaning 
6-dB  postsignal  processing  or  even  less)  and  are  then  applied  to 
five  synthetic  benchmark  multistatic  active  sonar  scenarios  fea¬ 
turing  multiple  targets,  multiple  sources,  and  multiple  receivers. 
Both  methods  end  up  performing  well  in  situations  where  there 
is  a  single  target  or  widely  spaced  targets.  However,  ML-PMHT 
has  an  inherent  advantage  over  ML-PDA  in  that  its  likelihood 
ratio  (LR)  has  a  simple  multitarget  formulation,  which  allows  it  to 
be  implemented  as  a  true  multitarget  tracker.  This  formulation, 
presented  here  for  the  first  time,  gives  ML-PMHT  superior  per¬ 
formance  for  instances  where  multiple  targets  are  closely  spaced 
with  similar  motion  dynamics. 

Index  Terms — Bistatic,  expectation  maximization  (EM),  low  ob¬ 
servable,  maximum  likelihood,  maximum-likelihood  probabilistic 
data  association  (ML-PDA),  maximum-likelihood  probabilistic 
multihypothesis  (ML-PMHT),  multistatic,  multitarget,  multi¬ 
target  ML-PMHT,  sonar,  tracking. 


I.  Introduction 

THE  maximum-likelihood  probabilistic  data  association 
(ML-PDA)  tracker  and  the  maximum-likelihood  proba¬ 
bilistic  multihypothesis  (ML-PMHT)  tracker  are  both  simple, 
straightforward  algorithms  that  can  be  used  in  an  active  mul¬ 
tistatic  sonar  framework.  With  some  basic  assumptions  about 
a  target  (or  targets)  as  well  as  the  environment,  likelihood 
ratios  (LRs)  can  be  developed  for  both  algorithms  and  then 
optimized.  The  main  difference  between  the  two  algorithms 
is  in  the  measurement-to-target  assignment  model;  ML-PDA 
assumes  that  at  most  one  measurement  per  scan  can  originate 
from  a  target,  while  ML-PMHT  allows  for  any  number  of 
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measurements  to  have  originated  from  a  target.  While  this 
assumption  may  reduce  the  appeal  of  ML-PMHT  to  some,  the 
resulting  algorithm  does  offer  advantages  in  both  its  implemen¬ 
tation  (especially  fine-scale  optimization)  and  in  terms  of  its 
multitarget  formulation. 

The  main  contributions  of  this  work  are  as  follows.  Lirst,  we 
theoretically  develop  and  implement  ML-PMHT  as  a  true  mul¬ 
titarget  tracker.  As  part  of  this,  we  find  an  expression  for  the 
Cramer-Rao  lower  bound  (CRLB)  for  ML-PMHT  and  show 
that  ML-PMHT  seems  to  be  a  statistically  efficient  estimator. 
As  a  result,  this  CRLB  can  be  used  as  a  covariance  estimate 
for  ML-PMHT.  Moving  on  to  Monte  Carlo  testing,  we  provide 
an  example  of  the  capability  of  the  algorithms  as  very  low  ob¬ 
servable  (VLO)  trackers.  [In  this  work,  we  use  “low  observ¬ 
able”  to  refer  to  targets  with  a  signal-to-noise  ratio  (SNR)  be¬ 
tween  6  and  12  dB  and  “very  low  observable”  to  refer  to  tar¬ 
gets  with  an  SNR  of  6  dB  or  less.]  We  then  demonstrate  that 
while  ML-PMHT  and  ML-PDA  have  identical  performance  in 
single-target  cases,  ML-PMHT  is  better  than  ML-PDA  in  mul¬ 
titarget  cases  or  when  more  than  one  measurement  per  scan  is 
produced  by  the  target. 

The  paper  has  the  following  structure.  Section  II  briefly  re¬ 
views  the  ML-PDA  and  ML-PMHT  (single-target)  log-like¬ 
lihood  ratios  (a  much  more  extensive  discussion  of  this  can 
be  found  in  [1]  and  [2]).  Section  III  develops  the  novel  multi¬ 
target  ML-PMHT  (including  the  CRLB  for  ML-PMHT).  Next, 
Section  IV  describes  in  detail  the  simulator  used  to  create  data 
for  Monte  Carlo  runs  on  the  multistatic  sonar  scenarios  con¬ 
sidered;  the  five  different  scenarios  used  in  this  work  are  pre¬ 
sented  here  as  well.  Linally,  Section  V  presents  the  results  of  the 
Monte  Carlo  testing.  We  note  that  portions  of  this  paper  have  ap¬ 
peared  in  [1]  and  [2];  this  paper  is  an  expansion  of  these  works. 
These  previous  works  provided  preliminary  comparisons  be¬ 
tween  ML-PDA  and  ML-PMHT  when  they  were  both  imple¬ 
mented  in  a  sequential  single-target  mode. 

II.  Single-Target  ML-PDA  Versus  ML-PMHT 

This  section  briefly  describes  the  theory  of  ML-PDA  and 
ML-PMHT  (single-target)  log-likelihood  ratios  (LLRs).  It  then 
shows  how,  under  certain  conditions,  the  ML-PDA  LLR  and 
the  ML-PMHT  LLR  should  converge  to  the  same  expression. 
If  such  conditions  are  met,  there  should  be  no  performance  dif¬ 
ference  between  the  two  algorithms. 
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A.  Likelihood  Ratios 

The  ML-PDA  concept  was  initially  developed  in  [3].  Subse¬ 
quently,  in  [4]-[6],  it  was  applied  in  a  multistatic  active  appli¬ 
cation,  which  is  how  we  currently  employ  it.  The  assumptions 
used  to  develop  the  ML-PDA  LR  are  [7],  [8]  as  follows. 

•  A  single  target  is  present  in  each  frame  with  known  de¬ 
tection  probability  Pd.  Detections  are  independent  across 
frames. 

•  There  are  zero  or  one  measurements  per  frame  from  the 
target. 

•  The  kinematics  of  the  target  are  deterministic.  The  motion 
is  usually  parameterized  as  a  straight  line,  although  any 
other  parameterizations  (e.g.,  ballistic  motion)  can  be  used. 

•  False  detections  (clutter)  are  uniformly  distributed  in  the 
search  volume  and  their  number  is  Poisson  distributed  with 
known  density. 

•  Amplitudes  of  target  and  false  detections  are  Rayleigh  dis¬ 
tributed.  The  parameter  of  each  Rayleigh  distribution  is 
known  (although  the  SNR  may  be  tracked  [9]  in  the  case 
that  it  is  not  known;  then  Pd  becomes  time  varying,  which 
is  easy  to  accommodate). 

•  Target  measurements  are  corrupted  with  zero-mean 
Gaussian  noise  with  known  variance. 

•  Measurements  at  different  times,  conditioned  on  the  pa¬ 
rameterized  state,  are  independent. 

The  assumptions  for  ML-PMHT  are  almost  all  the  same  as 
these,  with  one  significant  difference.  Instead  of  allowing  only 
zero  or  one  measurement  per  scan  to  originate  from  the  target, 
the  ML-PMHT  target-assignment  model  relaxes  this  constraint 
and  allows  any  number  of  measurements  to  originate  from  the 
target.  The  development  of  the  ML-PDA  and  ML-PMHT  LRs 
has  been  presented  in  detail  in  [2];  we  just  summarize  the  results 
here. 

The  ML-PDA  LLR  for  a  batch  of  data  is 


Nw 


i=  1 


H^z)  =  X]ln  { 1  ~  Pd  +  y  lx]PiW  )  •  C1) 

3  =  1 


In  this  equation,  Pd  is  the  target  probability  of  detection  in  a 
scan,  A  is  the  spatial  clutter  density,  Nw  is  the  number  of  scans, 
and  nrii  is  the  number  of  measurements  in  the  it h  scan.  Next, 
z j(i)  is  the  jth  measurement  in  the  it h  scan,  p[zj(i)  |  x]  is  a 
target-centered  Gaussian,  and  pj  (z)  is  the  amplitude  LR.  Finally, 
x  is  the  target  state,  and  Z  is  the  (entire)  set  of  measurements  in 
the  batch. 

The  ML-PMHT  LLR  is  given  by 

Nw  mi 

A'(x,  z)  =  +  ’’WpM*)  I  XW(*)}  (2) 

i=  1  j—  1 


where  7ro  is  the  prior  probability  that  a  measurement  is  from 
clutter,  7Ti  is  the  prior  probability  that  a  measurement  is  from 
the  target,  and  V  is  the  search  volume  in  the  measurement  space. 
Note  that  the  ML-PDA  LR  is  a  product  of  a  large  sum  (based 
on  the  total  probability  theorem  with  respect  to  all  the  joint  as¬ 
sociation  events)  while  the  ML-PMHT  LR  is  a  double  product 
of  the  sum  of  two  terms  (based  on  the  total  probability  theorem 
for  each  measurement  separately,  since  each  could  be  target  or 


clutter  originated,  and  these  events  are  assumed  independent 
across  measurements).  This  difference  is  the  key  in  1)  allowing 
the  computation  of  the  ML-PMHT  CRLB  in  a  much  easier 
manner  than  for  the  ML-PDA  CRLB;  and  2)  allowing  the  ex¬ 
tension  of  the  ML-PMHT  for  multiple  targets. 


B.  Relationship  Between  ML-PDA  and  ML-PMHT 

At  first  glance,  the  LLRs  for  ML-PDA  and  ML-PMHT 
appear  very  dissimilar.  However,  the  terms  in  the  LLR  for 
ML-PDA  are  actually  just  a  subset  of  the  terms  in  the  LR  for 
ML-PMHT.  Consider  the  ML-PMHT  likelihood  (not  log-like¬ 
lihood)  ratio  for  a  single  scan.  This  is  written  as 

!*j  =  jflUo  +  t’WpMO  I  x]Pj(*)}-  (3) 


For  any  arbitrary  number  of  measurements  rm,  this  equation 
multiplied  out  yields 


Pi[Z{i)  lx] 
Po[Z{i)  |x] 


3  =  1 

rri  i  —  1  m 

+ yt*  2'n\v2  53  53  + 

.71  -1  .72  >.7  i 

771  j 

+  7rm,ym;  j-J  yr 
3  =  1 


(4) 


where  the  terms  Afj  represent  p[zj(i)  |  x],  the  target-centered 
Gaussian  terms.  It  should  also  be  noted  that  the  amplitude  dis¬ 
tribution  ratios  have  been  left  out  of  the  above  equation;  in¬ 
cluding  them  would  not  affect  the  following  discussion.  Equa¬ 
tion  (4)  is  the  ML-PMHT  LR  multiplied  out  for  a  single  scan. 
It  allows  for  more  than  one  measurement  to  be  assigned  to  a 
target;  the  terms  with  products  of  target-centered  Gaussians  ac¬ 
count  for  this  case.  In  the  case  of  ML-PDA,  the  target  assign¬ 
ment  model  allows  only  zero  or  one  measurement  to  be  assigned 
to  the  target.  The  first  term  in  (4)  accounts  for  the  case  of  zero 
measurements  being  assigned  to  the  target,  and  the  second  term 
in  the  equation  accounts  for  the  case  of  one  measurement  being 
assigned  to  the  target.  The  other  terms  account  for  two  or  more 
measurements  being  assigned  to  the  target.  If  the  ML-PDA 
measurement  model  holds  (or  is  at  least  very  close  to  holding), 
then  the  product  of  the  probabilities  of  two  (or  more)  measure¬ 
ments  being  from  a  target  is  essentially  zero  and  the  terms  with 
a  product  of  two  or  more  target-centered  Gaussians  in  the  above 
expression  go  to  zero.  As  a  result,  just  the  first  and  second  terms 
are  left.  In  this  case,  we  can  claim  the  following  relationships. 
First 


1  -  Pd  =  7T 


771  i 

0  ' 


(5) 


Then 


7->  m  i  —  1 

Pd  =  mi7T0  7T1 


(6) 


and,  since  E[rrii\  =  XV,  one  can  take  (approximately) 

=  L.  (7) 

If  these  relationships  hold,  then  the  ML-PDA  LR  and  the 
ML-PMHT  LR  are  equivalent. 
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III.  Multitarget  ML-PMHT 

This  section  develops  the  multitarget  formulation  of 
ML-PMHT.  It  first  discusses  how  it  is  necessary  to  im¬ 
plement  ML-PDA  for  multiple  targets  in  somewhat  of  an  ad 
hoc  manner.  Next,  it  presents  the  natural,  multitarget  version 
of  the  ML-PMHT  LLR  along  with  an  elegant  method  of 
optimizing  it.  Then,  it  develops  a  covariance  estimate  of  the 
ML-PMHT  state  vector  via  the  CRLB.  Finally,  it  describes 
how  to  integrate  these  two  concepts  into  a  true  multitarget 
tracking  framework  for  ML-PMHT. 

A.  ML-PDA  Multitarget  LLR 

It  is  difficult  to  extend  ML-PDA  to  multiple  targets  [10]. 
While  it  is  technically  possible  to  write  the  multitarget 
ML-PDA  LLR,  to  take  into  account  all  the  joint  association 
events  the  number  of  terms  increases  rapidly  with  the  number 
of  targets,  and  the  expression  becomes  practically  intractable 
for  any  more  than  just  a  few  targets.  As  a  result,  to  handle 
multiple  targets,  the  ML-PDA  tracking  framework  treats  sce¬ 
narios  with  more  than  one  target  as  a  sequence  of  single-target 
problems.  For  a  batch  of  data,  it  optimizes  the  single-target 
LLR  (1),  and  if  this  value  exceeds  a  certain  threshold,  a  target 
is  declared.  (This  threshold  is  determined  by  fitting  the  local 
maxima  formed  by  the  clutter  to  an  extreme  value  distribution 
[11].)  Next,  the  measurement  that  has  the  highest  association 
probability  with  that  solution  is  excised  from  each  scan,  and 
the  sequence  is  repeated  for  the  next  target.  This  method  is 
not  elegant,  but  it  has  been  shown  to  work  reasonably  well  for 
multitarget  scenarios  in  [12]  and  [13]. 

B.  ML-PMHT  Multitarget  LLR 

In  contrast  to  ML-PDA,  the  ML-PMHT  LLR  is  very  easily 
extended  to  a  multitarget  framework.  For  K  targets  with  state 
vectors  xi, . . . ,  x^,  the  multitarget  LLR  is  expressed  as 

Nw  rrii  f  K  \ 

A'(*.Z)  =  £I>  7T0  +  V  ^  ^  |  X/C]pj^.(i)  /  (8) 

i=  1  j  =  1  L  k=  1  J 

where  itk  is  the  probability  that  a  given  measurement  is  from 
the  kth  target  and  7Tk  =  1- 

A  benefit  to  using  the  ML-PMHT  LLR  formulation  is  that 
it  is  possible  to  optimize  it  with  a  closed-form  expression 
using  expectation-maximization  (EM)  [14].  The  process  of 
optimizing  the  ML-PMHT  LLR  for  the  single-target  case  was 
shown  in  [1];  here,  we  extend  the  optimization  to  the  multitarget 
case.  The  multitarget  version  of  the  cost  function  is 

K  Nw  mz 

Ax,  z)  =  EEE  {M*)  -  HfcxfRUM*)  -  H*x] 

k= 1  i= 1  j= 1 

+  ln(|27rRij|)}  wjk(i).  (9) 

To  briefly  review,  in  the  single-target  case,  the  ML-PMHT  LLR 
can  be  optimized  with  EM  as  long  as  there  exists  a  linear  rela¬ 
tionship  between  the  predicted  measurement  z  and  the  state  x, 
where 


x  =  [xq  X  Uo  y]T 


and  the  measurement  matrix  is 


H  = 


1 

0 


t  0 
0  1 


0 

t 


(11) 


(For  this  to  work,  the  measurements  must  obviously  be  in  Carte¬ 
sian  space.  While  multistatic  measurements  and  covariances  are 
typically  in  bearing/time-delay  space,  they  can  be  converted  to 
Cartesian  space  following  the  work  of  [15]  and  [16].) 

In  the  multitarget  case  for  ML-PMHT,  the  state  vector  x  be¬ 
comes  a  4 K  x  1  vector.  [It  is  the  vector  in  (10)  stacked  K  times.] 
The  association  probability  Wjk(i)  is 


wjk{i) 


7Tfcp[Zj(i)  I  Xk]pjk(i) 

^o/V  +  YL=i  I 


Matrix  in  (9)  is 


H 


0 

0 

0 

0 


H 


0 

0 

0 

0 


4(fc  — 1)  columns  4 (K  —  k)  columns 


(12) 


(13) 


This  matrix  is  for  the  kth  target;  the  (inner)  H  is  that  given 
by  (1 1).  There  are  4(fc  —  1)  vectors  of  zeros  to  the  left  of  H  and 
4 (K  —  k)  vectors  of  zeros  to  the  right  of  H.  With  these  defi¬ 
nitions,  the  solution  to  the  problem  is  a  simple  vector  quadratic 
minimization,  and  is  expressed  as 


x  = 


i  y  w  m  i 

EEE 


K  Nw  mi 

xEEE  WjfcWHfeRAz^i). 

k= 1  i=l  j=  1 


(14) 


C.  ML-PMHT  CRLB 

Here,  the  CRLB  for  ML-PMHT  is  developed,  and  it  is  shown 
that  the  ML-PMHT  estimator  is  statistically  efficient,  i.e.,  the 
errors  on  the  tracker  state  estimates  meet  the  CRLB.  Addition¬ 
ally,  it  turns  out  that  the  CRLB  for  ML-PMHT  can  be  com¬ 
puted  in  real  time;  previous  work  on  ML-PDA  [8]  showed  that 
the  CRLB  for  ML-PDA  must  be  computed  offline  with  Monte 
Carlo  integration.  Put  together,  these  factors  lead  to  a  valuable 
result,  namely,  the  CRLB  calculation  can  be  easily  incorporated 
into  the  ML-PMHT  tracking  framework  to  represent  the  tracker 
error. 

The  full  derivation  for  the  ML-PMHT  Fisher  information 
matrix  (FIM)  and  CRLB  is  presented  in  the  Appendix;  we  sum¬ 
marize  the  key  points  to  show  why  it  is  possible  to  compute 
the  CRLB  in  real  time  for  ML-PMHT  (again,  as  opposed  to 
ML-PDA).  First,  consider  a  window  of  ML-PMHT  data  with 
Nw  scans.  Since  all  measurements  from  scan  to  scan  are  as¬ 
sumed  to  be  independent,  the  FIM  J  of  the  total  window  can 
just  be  written  as  the  sum  of  the  individual  scans 

Nw 

J  =  5>- 


(10) 


(15) 
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After  some  work  (again,  presented  in  the  Appendix),  a  key  point 
arises.  The  expression  for  J i  is  given  as  (16)  at  the  bottom 
of  the  page.  Here,  0  is  the  state-to-measurement  conversion, 
and  is  the  Jacobian  of  (j).  Next,  Rj  is  the  measurement  co- 
variance  for  the  jth  measurement  in  the  scan.  [Since  (16)  only 
deals  with  a  single  scan,  the  z-subscripts  are  dropped.]  Finally, 
dZ  =  Yli-i  dzi  and  da.  =  YYi=i  dai.  All  p(zi  |  x)  terms  for 
l  ^  j  are  separable  and  integrate  to  one,  and  the  term  for  l  =  j 
is  canceled  by  one  of  the  terms  in  the  denominator.  All  that  is 
left  within  the  integral  is  the  single  measurement  z7 .  This  re¬ 
duces  the  calculation  to  a  dimension  equal  to  the  number  of  di¬ 
mensions  in  the  measurement,  plus  one  for  the  amplitude.  This 
makes  it  possible  for  ML-PMHT  to  compute  the  FIM  (and,  thus, 
the  CRLB)  in  real  time  as  part  of  the  tracking  framework.  In  con¬ 
trast,  for  the  ML-PDA  FIM  calculation,  a  summation  over  all 
mi  measurements  in  the  scan  remains  in  the  denominator  [8]. 
Thus,  the  dimensionality  of  the  calculation  in  this  case  is  ap¬ 
proximately  mi  times  the  number  of  measurement  dimensions. 
As  a  result,  the  ML-PDA  FIM  must  be  calculated  offline  with 
Monte  Carlo  integration. 

After  some  work,  the  final  expression  for  the  ML-PMHT 
FIM  is 


j*  =  d£Eg 

3  = 1 


I, 


frlPlKQ]  ^  tT 

v  noPoM  +  KlPl(aj)  c-h£f£< 


G,  _ 

'GJ 


d£jda 


V 


V\  2?rRi 


(17) 


where  is  a  dummy  variable  of  integration  with  dimension¬ 
ality  equal  to  the  measurement  dimensionality,  and  is  the 
Cholesky  decomposition  of  R J1 .  This  is  inserted  into  (15),  and 
then  the  CRLB  is  just  taken  as  the  inverse  of  the  FIM. 

At  this  point,  the  efficiency  of  the  ML-PMHT  estimator  was 
checked  by  using  this  FIM  to  calculate  the  normalized  estima¬ 
tion  error  squared  (NEES),  which  is  given  by 

NEES  =  (x  -  xtruth)TJ(x  -  Xtruth)-  (18) 


Scenario  1  (described  below)  was  run  500  times,  and  for  each 
run,  the  median  NEES  value  for  a  valid  track  was  saved  (there 
were  multiple  NEES  values  for  a  run  due  to  the  sliding  batch  na¬ 
ture  of  the  tracker  implementation).  The  average  of  these  NEES 
values  came  out  to  4.15;  the  expected  NEES  for  a  four-param¬ 
eter  state  vector  is  4,  with  a  95%  confidence  interval  of  [3.81, 
4.19].  The  average  NEES  value  falls  within  this  confidence  in¬ 
terval,  so  the  conclusion  is  that  ML-PMHT  is  a  statistically  ef¬ 
ficient  estimator.  A  subset  of  20  of  these  runs  is  shown  in  Fig.  1 
with  the  95%  position  uncertainty  ellipses  at  the  beginning  and 


x  (arbitrary  units) 

Fig.  1.  CRLB  95%  uncertainty  regions  for  initial  and  final  positions. 


Bistatic  Angle  (degrees) 

Fig.  2.  TS  as  a  function  of  bistatic  angle. 


end  of  the  runs  plotted.  As  expected,  almost  all  of  the  solu¬ 
tion  points  fall  within  their  respective  ellipses.  The  CRLB  for 
ML-PMHT  (as  opposed  to  ML-PDA)  can  be  computed  in  real 
time  and  then  be  used  to  accurately  represent  the  error  on  the 
state  estimate  from  the  tracker. 

D.  ML-PMHT  Multitarget  Implementation 

The  ML-PMHT  multitarget  tracking  framework  was  im¬ 
plemented  in  the  track  update  sequence  by  testing  all  existing 
tracks  for  “closeness.”  Any  tracks  that  were  determined  to  be 
close  were  grouped  together  and  optimized  with  the  multitarget 
ML-PMHT  LLR  formulation  in  (8).  Grouping  tracks  involved 
estimating  the  state  covariance  of  all  existing  tracks.  We  let  the 
CRLB  developed  above  represent  the  state  covariance  Cn  for 
the  nth  target,  and  with  this,  a  x2  test  statistic  [17]  is  evaluated 


=  E  /,  D l - - D*  II  !>(*>  I 


3  = 1 


(16) 


i=i 
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Fig.  3.  All  measurements  from  a  window  of  data  (11  update  periods)  for  sce¬ 
nario  1  with  Rayleigh  clutter. 
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Fig.  4.  ML-PDA  in-track  percentage  versus  expected  target  SNR  for  scenario 
1  with  Rayleigh  clutter. 


for  closeness  between  all  possible  pairs  of  active  tracks  in  the 
form 


Smn  =  Ax^n(Cm  +  Cn)  1  Axmn  (19) 

where  Axmn  =  xm  —  xn  (the  difference  of  the  state  vector 
estimates  between  the  rath  and  nth  track).  All  track  pairs  with 
a  statistic  Smn  less  than  a  given  threshold  are  grouped  together. 
More  than  two  tracks  could  belong  to  a  group  by  the  use  of  this 
association;  to  join  a  group,  an  ungrouped  track  only  had  to  meet 
the  test  described  in  (19)  with  one  of  the  tracks  already  in  the 
group. 


IV.  Simulator  and  Simulations 

This  section  describes  in  detail  the  simulator  that  was  used  to 
create  the  data  for  the  Monte  Carlo  runs  and  then  presents  the 
five  scenarios  used  to  compare  ML-PDA  and  ML-PMHT. 

A.  Simulator  Description 

The  Monte  Carlo  runs  were  performed  using  a  multistatic 
simulator  developed  at  the  University  of  Connecticut  (Storrs, 
CT,  USA)  [18].  The  simulator  takes  as  input  the  source  and  re¬ 
ceiver  positions  (all  were  stationary  in  this  work),  as  well  as  the 
number  of  targets  and  target  geometries;  this  information  is  all 
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Fig.  5.  Scenario  1  (baseline)  and  ML-PDA  sample  estimated  track  in  one  run. 


x  (arbitrary  units) 

Fig.  6.  Scenario  1  (baseline)  and  ML-PMHT  sample  estimated  track  in  one 
run. 


shown  in  Figs.  2-14.  Each  transmitter  generates  a  ping  every 
60  s,  and  for  each  source-receiver  pair  (i.e.,  a  scan),  the  sim¬ 
ulator  calculates  the  number,  position,  amplitude,  and  Doppler 
of  the  clutter  points,  and  the  existence,  position,  amplitude,  and 
Doppler  of  the  target  returns.  (All  pings  in  the  simulation  were 
treated  as  being  Doppler  sensitive.) 

1)  Clutter  Measurement  Generation:  First,  to  generate  the 
clutter  points,  the  simulator  calculates  the  number  of  resolution 
cells  that  should  be  present.  The  number  of  resolution  cells  is 
treated  as  a  function  of  craz  and  at ,  the  azimuthal  and  time- 
delay  errors  simulated  for  the  system.  The  number  of  azimuthal 
resolution  cells  is  approximated  as 


iVa: 


360 

BW 


(20) 
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Fig.  7.  Scenario  2  (three  close  targets  with  similar  dynamics)  and  ML-PDA 
sample  estimated  tracks  in  one  run. 
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Fig.  9.  Scenario  3  (two  close  targets  with  different  dynamics)  and  ML-PDA 
sample  estimated  tracks  in  one  run. 


x  (arbitrary  units) 

Fig.  8.  Scenario  2  (three  close  targets  with  similar  dynamics)  and  ML-PMHT 
sample  estimated  tracks  in  one  run. 


where  BW  is  the  beamwidth  of  a  receive  beam  in  degrees.  The 
simulator  does  not  actually  form  any  beams;  to  get  a  value  of 
BW,  it  makes  the  assumption  that  the  clutter  is  distributed  uni¬ 
formly  in  a  beam,  which  gives  craz  =  B W / \/l2 .  This,  in  turn, 
leads  to  a  final  expression  for  the  number  of  azimuthal  resolu¬ 
tion  cells 


iVaz  = 


360 


a/12  o 


(21) 


(Here,  signifies  round  toward  the  next  lower  integer.)  For 
all  of  the  scenarios  in  this  work,  the  azimuthal  uncertainty  was 
<raz  =  5°  (this  and  all  other  simulated  error  values  are  provided 
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Fig.  10.  Scenario  3  (two  close  targets  with  different  dynamics)  and  ML-PMHT 
sample  estimated  tracks  in  one  run. 


together  in  Table  II).  Using  similar  logic,  the  number  of  delay 
(range)  resolution  cells  is  calculated  with  the  expression 


^delay  — 


^ipt/2  ^bt 

VU(Tt 


(22) 


Here,  tipt  is  inter-ping  time  (60  s  throughout  this  work),  and  t^t 
is  blanking  time.  This  value  is  determined  by  the  amount  of  time 
sound  takes  to  travel  from  the  source  to  the  receiver.  Finally,  the 
number  of  resolution  cells  is  simply  given  by  the  product  of  7Vaz 

and  N delay  • 

With  the  number  of  resolution  cells  determined,  the  clutter 
is  generated.  For  each  scan,  to  simulate  the  amplitudes  of  the 
clutter  distribution,  a  number  of  random  variables  equal  to  the 
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x  (arbitrary  units) 

Fig.  11.  Scenario  4  (large  number  of  targets)  and  ML-PDA  sample  estimated 
tracks  in  one  run. 


Fig.  12.  Scenario  4  (large  number  of  targets)  and  ML-PMHT  sample  estimated 
tracks  in  one  run. 


Fig.  13.  Scenario  5  (switching  targets)  and  ML-PDA  sample  estimated  tracks 
in  one  run  (switch  occurred). 


x  (arbitrary  units) 

Fig.  14.  Scenario  5  (switching  targets)  and  ML-PMHT  sample  estimated 
tracks  in  one  run,  zoomed-in  view  (no  switching  occurred). 

Random  variables  from  the  distribution  shown  in  (23)  are  gen¬ 
erated  via  [20] 


number  of  resolution  cells  is  generated.  For  the  first  five  sce¬ 
narios,  the  amplitude  of  this  clutter  distribution  was  given  a 
K -distribution.  The  probability  distribution  for  the  intensity  of 
the  clutter  (i.e.,  the  amplitude  squared)  is  given  by  [19],  [20] 


f(p)  = 


2 

ATK) 


(afc-l)/2 


Kak-\ 


(23) 


where  p  >  0.  Here,  is  a  K -distribution  parameter  (ak  = 
0.5  throughout  this  work),  Xk  =  1  /ak,  and  K„  is  the  Basset 
function  (a  modified  Bessel  function  of  the  second  kind)  [21]. 


K  =  -TXk\nU  (24) 

where  U  is  a  realization  of  a  uniform  random  variable  over  the 
interval  [0,  1],  and  T  is  a  realization  of  a  Gamma  random  vari¬ 
able  with  shape  parameter  a  =  ak  and  a  scale  parameter  f3  =  1. 

For  clarity,  we  note  that  the  simulator  generates  “amplitude” 
LLRs  by  operating  in  intensity  (power)  space,  i.e.,  amplitude 
squared.  Such  LLRs  can  be  computed  in  either  amplitude  or 
intensity  space;  as  long  as  the  clutter  and  target  distributions  are 
consistent  in  terms  of  the  units  used,  the  LLR  values  will  come 
out  the  same. 
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TABLE  I 

Detector  Thresholds 


Scenario 

Threshold  (dB) 

1-4 

12 

5 

14 

1  (Rayleigh  clutter) 

2 

TABLE  II 

Simulation  Errors 


&az 

(deg) 

O'* 

(S) 

O’Doppler  target 

(units/s) 

Doppler  clutter 

(units/s) 

5 

0.1 

0.5 

2 

For  the  final  scenario,  the  clutter  amplitude  was  given  a 
Rayleigh  distribution  following  the  work  in  [8],  so  the  intensity 
of  the  clutter  followed  an  exponential  distribution.  The  clutter 
intensity  in  this  case  was  simulated  by  generating  exponential 
random  variables  with  a  mean  of  one. 

Once  the  clutter  realizations  were  created,  they  were  thresh- 
olded.  Amplitudes  above  the  threshold  were  kept;  amplitudes 
below  the  threshold  were  discarded.  This  threshold  is  what  de¬ 
termined  the  clutter  density;  it  was  set  so  that  the  clutter  density 
for  the  average  source-receiver  pair  in  the  first  five  scenarios 
was  approximately  1  x  10“ 9  clutter  points  per  unit  volume;  for 
the  final  scenario  with  Rayleigh  clutter,  the  threshold  was  set 
so  that  the  clutter  density  was  approximately  1  x  10  _  7  clutter 
points  per  unit  volume.  The  amplitude  detector  threshold  values 
(in  decibels)  are  given  in  Table  I.  After  this,  all  remaining  clutter 
points  were  given  delay  time  and  azimuth  measurements  by 
sampling  from  a  uniform  distribution  over  the  proper  time  and 
azimuthal  interval.  Finally,  all  generated  clutter  measurements 
were  given  a  Doppler  component.  The  assumption  was  made 
that  the  clutter  would  be  stationary;  thus,  the  Doppler  measure¬ 
ments  (simulated  in  the  form  of  bistatic  range  rate)  were  gen¬ 
erated  from  a  zero-mean  Gaussian  with  <rDoppler  ciutter  =  2 
units/s. 

2)  Target  Measurement  Generation:  For  every  update,  the 
target  and  source  positions  were  updated  according  to  user-de¬ 
fined  course  and  speed  segments.  The  expected  target-reflected 
power  at  the  receiver  was  then  calculated  (in  decibels)  with 

Prec  =  SLiooo  -  101og10(f^sT/1000) 

-101og10(fWl000)  +  10  log10(TS).  (25) 

Here,  SLiooo  is  the  signal  power  1000  distance  units  from  the 
target  (all  distance  units  in  this  simulation  are  arbitrary),  RT g 
is  the  distance  from  the  target  to  the  source,  Rtr  is  the  dis¬ 
tance  from  the  target  to  the  receiver,  and  TS  is  the  bistatic  target 
strength.  This  bistatic  target  strength  is  a  function  of  the  bistatic 
angle,  defined  in  [16].  For  this  simulation,  the  value  of  TS  as 
a  function  of  bistatic  angle  is  shown  in  Fig.  2.  [Note  that  with 
(25)  we  are  not  trying  to  implement  directly  the  physical  effects 
described  by  the  active  sonar  equation.  Instead,  we  are  trying 
to  implement  a  target  amplitude  model  that  has  a  known  de¬ 
pendency  on  range  and  bistatic  aspect  and  is  given  a  reference 
level  at  a  convenient,  fixed  distance  from  the  target.]  Finally,  the 
“measured”  target  amplitude  was  simulated  having  a  Rayleigh 


TABLE  III 

SLiooo  Values  as  a  Function  of  Scenario 


Scenario 

1 

2 

3 

4 

5 

1  (Rayleigh  clut) 

SLiooo  (dB) 

35 

40 

45 

40 

45 

27 

distribution,  with  expected  amplitude  given  by  (25).  The  mea¬ 
sured  target  amplitude  random  variable  is  simulated  (in  decibel 
space)  with 

-Ltgt  meas  —  Prec  +  101og10(— In  ([/)).  (26) 

Here,  U  is  a  uniform  random  variable  in  the  interval  [0,  1].  [It 
is  easy  to  show  that  if  Ptgt  meas  is  converted  to  intensity  space, 
(26)  produces  a  realization  of  an  exponential  random  variable 
with  mean  ioPrec/10 .]  If  target  return  power  Ptgt  meas  is  greater 
than  the  respective  threshold  for  its  scenario  shown  in  Table  I, 
then  the  (target-originated)  detect  is  kept  in  the  data  for  pro¬ 
cessing.  The  values  of  SLiooo  are  shown  in  Table  III. 

Finally,  the  time,  azimuth,  and  bistatic  range-rate  values  for 
all  target  measurements  are  calculated  and  then  corrupted  with 
noise.  Based  on  the  (deterministic)  positions  of  the  source,  the 
receiver,  and  the  target,  delay  time  t  and  azimuth  0  (from  the  re¬ 
ceiver  to  the  target)  are  determined  via  simple  geometry.  These 
values  of  t  and  0  are  corrupted  by  adding  realizations  of  zero- 
mean  Gaussian  random  variables  with  variances  of  and  craz 
to  them,  respectively.  The  bistatic  range  rate  is  calculated  fol¬ 
lowing  the  work  of  [16],  and  then  this  is  corrupted  with  a  zero- 
mean  Gaussian  random  variable  with  variance  of)oppler  target . 

3)  ML-PMHT  Target  Measurement  Generation  Model:  The 
majority  of  simulations  in  this  work  followed  the  ML-PDA 
target  measurement  generation  model,  that  is,  at  most  one  mea¬ 
surement  originated  from  the  target  in  a  given  scan  (this  is  im¬ 
plicit  in  the  simulator  discussion  above).  However,  the  simulator 
was  also  set  up  with  the  ability  to  generate  multiple  target  mea¬ 
surements  per  scan,  i.e.,  the  ML-PMHT  target  measurement 
generation  model.  This  was  done  by  assuming  Pd  is  known,  and 
then  writing  the  probability  mass  function  (PMF)  describing  the 
number  of  target  measurements  f(x.\t)  as 


I’  1  -Pd, 

x  =  0 

f(x,Xt)  =  < 

P,J  - 

x  >  0. 

(27) 

l  1  -  e  V 

x\ 

Let  N  =  E[x ],  the  expected  number  of  measurements  from 
the  PMF.  Calculating  the  expected  value  of  (27)  produces  the 
following  relationship: 

N(l-e~x')=Pd\t.  (28) 

Again,  the  value  of  Pd  is  assumed  known  (in  this  work,  it  was 
Pd  ~  0.8).  Then,  (28)  can  be  numerically  solved  for  \t  as  a 
function  of  N.  Results  of  this  are  presented  in  Table  IV.  To  pro¬ 
duce,  on  average,  N  target-generated  measurements  per  scan, 
the  simulator  first  determines  if  the  contact  power  is  greater  than 
the  threshold,  as  described  above.  If  the  target  power  does  ex¬ 
ceed  the  threshold,  the  simulator  draws  a  realization  from  a  trun¬ 
cated  Poisson  PMF  with  parameter  \t ,  and  creates  this  many 
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TABLE  IV 

Values  for  \t  as  a  Function  of  the  Expected 
Number  of  Target  Measurements 


N 

At 

1 

0.4642 

2 

2.2316 

3 

3.5628 

10 

7.5 

5 
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Fig.  15.  Scenario  1  (baseline)  and  ML-PDA  sample  estimated  track  with  high 
density  Rayleigh  clutter  in  one  run. 


(target)  measurements,  again,  in  the  same  manner  as  described 
above. 

B.  Scenarios 

Five  benchmark  multistatic  sonar  scenarios,  initially  devel¬ 
oped  in  [1],  were  used  to  measure  performance  differences 
between  ML-PDA  and  ML-PMHT  with  Monte  Carlo  testing. 
Each  scenario  was  designed  so  target  detections  from  a  single 
source-receiver  pair  would  be  present  approximately  80%  of 
the  time  in  a  given  scan,  and  clutter  was  set  at  a  level  that 
made  the  problem  as  difficult  as  possible  while  not  slowing 
down  run  times  to  the  point  of  precluding  Monte  Carlo  testing. 
The  various  scenario  parameters  are  listed  in  Table  V  (these 
parameters,  used  in  the  simulation,  were  also  matched  in  the 
actual  ML-PDA  and  ML-PMHT  tracking  code).  All  scenarios 
are  shown  (with  example  results  overlaid)  in  Figs.  3-16.  Each 
of  the  scenario  target  geometries  was  designed  for  a  specific 
test  purpose.  Then,  given  these  target  geometries,  the  sensors 
were  laid  out  in  an  effort  to  obtain  an  average  target  Pd  of  0.8 
in  any  given  scan.  (The  one  exception  to  this  occurs  in  scenario 
5  and  is  described  below.)  The  pinging  strategy  was  purposely 
kept  very  simple;  every  transmitter  had  a  (simulated)  ping 
every  60  s.  The  scenarios  are  described  as  follows. 

1)  Baseline  Scenario:  Scenario  1,  shown  in  Fig.  3,  featured 
a  single  target  moving  in  a  straight  line  past  a  source  and  a  re¬ 
ceiver  (the  source  was  a  receiver  as  well).  This  was  intentionally 
created  as  one  of  the  simplest  possible  multistatic  scenarios.  As 


Fig.  16.  Scenario  1  (baseline)  and  ML-PMHT  sample  estimated  track  with 
high  density  Rayleigh  clutter  in  one  run. 


TABLE  V 

ML-PMHT  AND  ML-PDA  PARAMETERS 


ML-PMHT 

ML-PDA 

7T0 

7T1 

V 

A 

Pd 

Seen.  1-4 

0.95 

0.05 

1.26  x  109 

3.7  x  10“9 

0.8 

Seen.  5 

0.95 

0.05 

1.26  x  109 

6.7  x  10-10 

0.8 

Seen.  1 

Rayleigh  clut. 

0.95 

0.05 

1.26  x  109 

2.0  x  10-7 

0.8 

is  shown  in  Section  II,  ML-PDA  and  ML-PMHT  should  con¬ 
verge  to  the  same  value  under  benign  conditions,  and  this  sce¬ 
nario  is  a  case  where  this  convergence  should  hold  true. 

2)  Close  Targets  With  Similar  Dynamics  Scenario:  Scenario 
2,  shown  in  Fig.  7,  features  three  targets  very  close  together 
(they  have  a  separation  of  only  500  distance  units)  with  similar 
motion  dynamics,  i.e.,  they  are  maneuvering  in  a  coordinated 
fashion.  This  makes  it  very  difficult  for  any  algorithm  to  distin¬ 
guish  and  separate  targets  from  each  other.  This  will  test  if  the 
multitarget  implementation  for  ML-PMHT  is  an  improvement 
over  the  sequential  single-target  implementation  of  ML-PDA. 

3)  Close  Targets  With  Different  Dynamics  Scenario:  Sce¬ 
nario  3,  shown  in  Fig.  9,  has  two  targets  closing  each  other  from 
the  east  and  west.  During  the  middle  of  the  scenario,  the  targets 
will  be  in  close  proximity  to  each  other,  but  they  will  have  dif¬ 
ferent  motion  dynamics;  at  the  time  they  are  close,  they  will  be 
moving  in  opposite  directions.  Again,  this  is  designed  to  test 
each  tracker’s  ability  to  follow  close  targets,  but  under  slightly 
easier  circumstances  than  in  scenario  2. 

4)  Large  Number  of  Targets  Scenario:  Scenario  4,  shown  in 
Fig.  11,  features  ten  very  low-speed  targets  and  three  relatively 
high-speed  targets.  This  scenario  was  designed  simply  to  mea¬ 
sure  each  algorithm’s  ability  to  track  a  (relatively)  large  number 
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of  targets,  especially  when  some  of  those  targets  are  low  speed 
and  exhibit  very  low  Doppler. 

5)  Switching  Targets  Scenario :  Scenario  5,  shown  in  Fig.  13, 
is  similar  to  scenario  2  in  that  it  has  two  targets  in  close  prox¬ 
imity  to  each  other  with  the  same  motion  dynamics.  However, 
the  targets  start  out  with  a  distinct  separation  and  then  close  on 
each  other.  At  the  point  where  the  targets  are  about  to  intersect, 
they  turn  and  parallel  each  other.  Such  motion,  with  the  targets 
approaching  each  other  (where  their  projected  motion  has  them 
crossing),  will  make  it  very  difficult  for  any  trackers  continu¬ 
ally  to  associate  measurements  with  the  correct  targets  and  not 
to  switch  targets.  Additionally,  this  scenario  was  given  a  rela¬ 
tively  high  number  of  receivers  (16),  and  due  to  the  geometry, 
most  source  receiver  pairs  will  have  a  Pd  of  much  less  than  0.8; 
many  scans  will  just  have  clutter  measurements.  This  is  effec¬ 
tively  raising  the  clutter  levels  that  the  tracker  will  see. 

These  scenarios  were  all  run  using  the  ML-PDA  target  mea¬ 
surement  generation  model,  where  the  target  produced  at  most 
one  measurement  per  scan.  Scenario  1  was  then  rerun  with  the 
ML-PMHT  target  measurement  generation  model,  where  more 
than  one  measurement  was  allowed  to  be  generated  by  the  target 
in  a  scan.  This  scenario  was  run  three  times  following  the  target 
measurement  generation  method  described  in  Section  IV- A3; 
the  expected  number  of  target  measurements  per  scan  was  one, 
two,  and  three,  respectively.  (Note  the  case  of  one  target  mea¬ 
surement  per  scan  is  not  the  ML-PDA  target  measurement  gen¬ 
eration  model.  In  this  case,  on  average,  the  target  generates  one 
measurement  per  scan,  but  there  could  be  more.  In  the  ML-PDA 
case,  the  number  of  target  measurements  per  scan  is  strictly  lim¬ 
ited  to  zero  or  one.) 

At  the  conclusion  of  each  scenario,  the  following  metrics 
were  evaluated:  target  in-track  percentage,  target  position  root 
mean  square  error  (RMSE),  track  fragmentation,  target  dupli¬ 
cate  tracks,  overall  number  of  false  tracks,  and  average  false 
track  length.  These  metrics  were  calculated  in  the  following 
manner.  First,  tracks  were  associated  with  either  a  target  or 
clutter.  A  track  was  associated  with  a  target  if  its  NEES  value 
was  less  than  10  for  scenarios  with  a  single  target  or  widely 
spaced  targets.  (In  most  situations,  as  is  discussed  above,  the  es¬ 
timators  were  statistically  efficient.)  For  scenarios  with  closely 
spaced  maneuvering  targets,  the  CRLB  tended  to  underestimate 
the  covariance  (the  estimators  were  not  statistically  efficient 
here),  so  the  NEES  was  not  as  reliable  in  associating  tracks  to 
targets.  Instead,  in  these  cases,  a  track  was  associated  with  a 
target  if  the  average  distance  between  track  points  and  their  as¬ 
sociated  truth  points  was  less  than  2000  distance  units.  If  a  track 
could  be  associated  with  more  than  one  target,  then  it  was  as¬ 
signed  to  whichever  one  was  closest.  Target  in-track  percentage 
was  then  calculated  by  taking  the  ratio  of  target  truth  points  that 
had  a  track  associated  with  them  to  the  total  number  of  target 
truth  points.  Target  duplicate  tracks  were  calculated  by  simply 
counting  the  number  of  targets  that  had  more  than  one  track  as¬ 
sociated  with  them  at  a  given  time.  The  RMSE  was  then  calcu¬ 
lated  for  all  tracks  that  were  associated  with  targets.  This  RMSE 
value  was  only  calculated  for  a  target  on  points  where  there  were 
both  an  ML-PDA  track  and  an  ML-PMHT  track  associated 
with  it.  This  was  done  so  as  not  to  penalize  one  algorithm  that 
was  barely  holding  track  (with  a  resultant  poor  RMSE)  while 


TABLE  VI 

In-Track  Percentage  Results 


ML-PDA 

ML-PMHT 

mean 

confidence  interval 
(95  percent) 

mean 

confidence  interval 
(95  percent) 

Scenario  1 

95.4 

[94.8,  96.0] 

95.4 

[94.8,  96.0] 

Scenario  2  Tgt  1 

12.0 

[7.5,  16.6] 

83.9 

[79.6,  88.2] 

Scenario  2  Tgt  2 

100 

[100,  100] 

87.5 

[83.7,  91.3] 

Scenario  2  Tgt  3 

1.9 

[0.0,  3.8] 

53.1 

[46.5,  59.7] 

Scenario  3  Tgt  1 

82.7 

[81.0,  84.4] 

82.5 

[81.0,  84.0] 

Scenario  3  Tgt  2 

76.1 

[74.2,  78.0] 

77.1 

[75.4,  78.9] 

Scenario  4  Tgt  7 

61.0 

[54.7,  67.4] 

91.1 

[89.3,  92.8] 

Scenario  4  Tgt  12 

52.5 

[45.9,  59.2] 

96.3 

[95.9,  96.8] 

Scenario  5  Tgt  1 

45.3 

[38.4,  52.2] 

63.1 

[56.5,  69.7] 

Scenario  5  Tgt  2 

48.4 

[41.5,  55.4] 

63.7 

[57.1,  70.2] 

Scenario  1  Rayleigh 

66.9 

[66.0,  67.8] 

67.4 

[66.6,  68.3] 

TABLE  VII 
RMSE  Results 


ML-PDA 

ML-PMHT 

mean 

confidence  interval 
(95  percent) 

mean 

confidence  interval 
(95  percent) 

Scenario  1 

255.8 

[229.7,  281.9] 

259.2 

[233.0,  285.4] 

Scenario  2  Tgt  1 

313.4 

[257.7,  369.0] 

281.4 

[267.7,  295.1] 

Scenario  2  Tgt  2 

247.0 

[241.9,  252.1] 

280.5 

[269.6,  291.5] 

Scenario  2  Tgt  3 

228.4 

[186.7,  270.0] 

212.4 

[202.3,  222.4] 

Scenario  3  Tgt  1 

1227.3 

[1190.2,  1264.3] 

946.6 

[904.4,  988.8] 

Scenario  3  Tgt  2 

799.7 

[754.6,  844.8] 

878.0 

[836.1,  920.0] 

Scenario  4  Tgt  7 

391.3 

[303.8,  478.8] 

321.6 

[274.5,  377.7] 

Scenario  4  Tgt  12 

151.6 

[124.0,  179.2] 

213.2 

[177.2,  249.3] 

Scenario  5  Tgt  1 

213.1 

[173.3,  252.8] 

244.7 

[207.1,  282.2] 

Scenario  5  Tgt  2 

272.5 

[211.3,  333.7] 

245.6 

[215.2,  276.1] 

Scenario  1  Rayleigh 

217.6 

[191.7,  243.6] 

182.4 

[166.4,  198.9] 

the  other  algorithm  was  not  tracking  at  all.  Track  fragmentation 
was  determined  by  counting  the  number  of  breaks  in  track  for 
a  target.  Finally,  the  overall  number  of  false  tracks  was  simply 
the  number  of  tracks  that  were  not  associated  with  a  target,  and 
the  mean  false  track  length  was  the  average  number  of  updates 
for  which  a  false  track  was  active. 

For  each  of  the  Monte  Carlo  sets,  200  runs  were  performed. 
Additionally,  the  confidence  intervals  shown  in  Tables  VI-XII 
are  at  the  95%  level. 

V.  ML-PDA  Versus  ML-PMHT  Performance 
Comparisons 

Both  ML-PDA  and  ML-PMHT  showed  themselves  to  be 
excellent  trackers  for  targets  with  low  amplitude  returns  (de¬ 
tailed  below).  First,  to  get  an  idea  of  the  level  of  difficulty  of  the 
problem,  all  measurements  from  a  batch  are  plotted  in  Fig.  3  for 
scenario  1  with  the  Rayleigh  clutter  distribution  (the  scenario 
with  the  highest  level  of  clutter). 

Using  this  same  scenario,  to  provide  a  comparison  to  other 
algorithms  with  known  performance,  in-track  percentage  as  a 
function  of  SNR  was  calculated.  This  is  shown  in  Fig.  4  for 
ML-PDA  (results  were  practically  identical  for  ML-PMHT). 
For  this  example  (in  contrast  to  the  rest  of  the  runs),  expected 
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TABLE  VIII 

Fragmentation  Results 


ML-PDA 

ML-PMHT 

mean 

confidence  interval 
(95  percent) 

mean 

confidence  interval 
(95  percent) 

Scenario  1 

0.01 

[0.00,  0.02] 

0.01 

[0.00,  0.02] 

Scenario  2  Tgt  1 

0.00 

[0.00,  0.00] 

0.00 

[0.00,  0.00] 

Scenario  2  Tgt  2 

0.00 

[0.00,  0.00] 

0.00 

[0.00,  0.00] 

Scenario  2  Tgt  3 

0.00 

[0.00,  0.00] 

0.00 

[0.00,  0.00] 

Scenario  3  Tgt  1 

0.03 

[0.01,  0.06] 

0.07 

[0.03,  0.10] 

Scenario  3  Tgt  2 

0.24 

[0.17,  0.30] 

0.06 

[0.02,  0.09] 

Scenario  4  Tgt  7 

0.00 

[0.00,  0.00] 

0.00 

[0.00,  0.00] 

Scenario  4  Tgt  12 

0.00 

[0.00,  0.00] 

0.00 

[0.00,  0.00] 

Scenario  5  Tgt  1 

0.15 

[0.07,  0.22] 

0.15 

[0.09,  0.22] 

Scenario  5  Tgt  2 

0.09 

[0.03,  0.15] 

0.15 

[0.08,  0.22] 

Scenario  1  Rayleigh 

0.07 

[0.03,  0.10] 

0.05 

[0.02,  0.08] 

TABLE  IX 

Duplicate  Track  Results 


ML-PDA 

ML-PMHT  | 

mean 

confidence  interval 
(95  percent) 

mean 

confidence  interval 
(95  percent) 

Scenario  1 

0.18 

[0.13,  0.24] 

0.02 

[0.00,  0.04] 

Scenario  2  Tgt  1 

0.00 

[0.00,  0.00] 

0.04 

[0.01,  0.08] 

Scenario  2  Tgt  2 

1.07 

[1.01,  1.14] 

0.48 

[0.39,  0.57] 

Scenario  2  Tgt  3 

0.00 

[0.00,  0.00] 

0.00 

[0.00,  0.00] 

Scenario  3  Tgt  1 

0.03 

[0.01,  0.06] 

0.07 

[0.03,  0.10] 

Scenario  3  Tgt  2 

0.24 

[0.17,  0.30] 

0.06 

[0.02,  0.09] 

Scenario  4  Tgt  7 

0.00 

[0.00,  0.00] 

0.09 

[0.05,  0.13] 

Scenario  4  Tgt  12 

0.09 

[0.03,  0.14] 

0.18 

[0.12,  0.24] 

Scenario  5  Tgt  1 

0.04 

[0.00,  0.08] 

0.13 

[0.07,  0.19] 

Scenario  5  Tgt  2 

0.02 

[0.00,  0.05] 

0.11 

[0.05,  0.16] 

Scenario  1  Rayleigh 

2.10 

[1.83,  2.37] 

0.00 

[0.00,  0.00] 

table  X 

False  Track  Results 


ML-PDA 

ML-PMHT  ] 

mean 

confidence  interval 
(95  percent) 

mean 

confidence  interval 
(95  percent) 

Scenario  1 

0.07 

[0.03,  0.11] 

0.07 

[0.03,  0.11] 

Scenario  2 

0.06 

[0.03,  0.09] 

0.15 

[0.15,  0.21] 

Scenario  3 

0.34 

[0.26,  0.42] 

0.26 

[0.19,  0.33] 

Scenario  4 

0.17 

[0.10,  0.23] 

0.14 

[0.09,  0.20] 

Scenario  5 

2.13 

[1.92,  2.34] 

2.65 

[2.44,  2.86] 

Scenario  1  Rayleigh 

0.08 

[0.02,  0.14] 

0.00 

[0.00,  0.00] 

target  SNR  was  fixed,  starting  at  4  dB  with  a  2-dB  measure¬ 
ment  threshold.  The  expected  target  SNR  was  then  moved  up 
in  1-dB  increments,  with  Monte  Carlo  simulations  being  run  at 
each  SNR  level.  Note  that  the  tracker  achieves  a  Pd  of  approx¬ 
imately  50%  at  a  target  SNR  of  only  4  dB,  and  a  Pd  of  about 
90%  at  6  dB.  By  a  target  SNR  of  8  dB,  the  tracker  essentially 
has  a  Pd  of  100%.  This  shows  that  ML-PDA  and  ML-PMHT 
are  effective  VLO  trackers. 


table  xi 

False  Track  Length  Results 


ML-PDA 

ML-PMHT 

mean 

confidence  interval 
(95  percent) 

mean 

confidence  interval 
(95  percent) 

Scenario  1 

3.14 

[2.59,  3.69] 

3.14 

[2.56,  3.73] 

Scenario  2 

2.75 

[2.31,  3.19] 

3.37 

[2.81,  3.92] 

Scenario  3 

4.34 

[3.18,  5.50] 

4.42 

[3.32,  5.62] 

Scenario  4 

6.75 

[3.72,  9.78] 

11.67 

[7.81,  15.52] 

Scenario  5 

8.98 

[7.81,  10.16] 

8.55 

[7.59,  9.51] 

Scenario  1  Rayleigh 

3.58 

[2.39,  4.78] 

0.00 

[0.00,  0.00] 

table  XII 

Number  of  Duplicate  Tracks  as  a  Function  of  N,  the  Expected 
Number  of  Target  Measurements  Per  Scan 


ML-PDA 

ML-PMHT 

N 

mean 

confidence  interval 
(95  percent) 

mean 

confidence  interval 
(95  percent) 

1 

1.32 

[1.21,  1.43] 

0.03 

[0.00,  0.05] 

2 

4.05 

[3.88,  4.23] 

0.04 

[0.01,  0.06] 

3 

6.12 

[5.92,  6.31] 

0.09 

[0.05,  0.13] 

After  this,  ML-PDA  and  ML-PMHT  were  run  on  all  five  sce¬ 
narios  in  their  respective  multitarget  tracking  modes.  ML-PDA 
operated  in  the  sequential  single-target  mode,  while  ML-PMHT 
was  run  in  its  true  multitarget  mode  for  tracks  that  were  close 
together. 

All  the  scenarios  except  one  featured  clutter  that  had  K- dis¬ 
tributed  amplitudes.  The  high-clutter  version  of  scenario  1  fea¬ 
tured  Rayleigh-distributed  clutter  amplitudes.  In  all  cases,  the 
feature  LLRs  calculated  by  ML-PDA  and  ML-PMHT  matched 
the  correct  clutter  distributions. 

For  scenarios  with  a  single  target  or  easily  distinguished  tar¬ 
gets,  ML-PDA  and  ML-PMHT  had  virtually  identical  perfor¬ 
mance  over  the  Monte  Carlo  runs.  As  is  shown  in  Section  II-B, 
the  ML-PMHT  LLR  should  be  very  close  to  the  ML-PDA  LLR, 
so  the  tracking  results  of  the  two  should  be  similar.  The  terms 
in  (4)  that  represent  more  than  one  measurement  in  a  scan  origi¬ 
nating  from  the  target  are  clearly  very  close  to  zero,  resulting  in 
the  ML-PDA  and  ML-PMHT  LLRs  being  nearly  identical  and 
thus  producing  nearly  identical  tracking  results. 

Individual  examples  taken  from  the  Monte  Carlo  runs  help  to 
further  illustrate  this.  For  scenario  1,  an  ML-PDA  result  shown 
in  Fig.  5  and  an  ML-PMHT  result  shown  in  Fig.  6  are  virtually 
indistinguishable  from  each  other.  For  the  same  scenario  with 
the  Rayleigh  clutter,  results  for  the  two  algorithms  (Fig.  15  for 
ML-PDA  and  Fig.  16  for  ML-PMHT)  are  also  very  similar. 
Finally,  for  scenario  3,  the  ML-PDA  example  in  Fig.  9  and  the 
ML-PMHT  example  in  Fig.  10  appear  very  much  alike. 

For  the  scenarios  with  targets  in  close  proximity  to  each 
other  with  similar  motion  dynamics  (scenarios  2  and  5),  or  the 
scenario  with  multiple  targets  that  were  close  enough  to  poten¬ 
tially  cause  a  tracker  to  switch  targets  (scenario  4),  the  results 
in  Tables  VI-XI  show  that  ML-PMHT  clearly  outperforms 
ML-PDA.  This  is  due  to  the  true  multitarget  implementation 
for  ML-PMHT  in  contrast  to  the  sequential  single-target 
tracking  logic  that  is  necessary  for  ML-PDA.  In  scenario  2, 
ML-PMHT  is  able  to  track  targets  1  and  2  in  excess  of  80%  of 
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•  High  SNR  measurement  •  Medium  SNR  measurement 
Low  SNR  measurement 


Target  Motion  (for  all  three  targets) - ^ 


Measurements  |  |  Measurements 

excised  by  track  1  I — I  excised  by  track  2 

- Track  1  - Track  2 


Fig.  17.  Example  of  multitrack  measurement  assignment  for  ML-PDA. 

the  time  and  track  target  3  nearly  50%  of  the  time.  In  contrast, 
ML-PDA  is  tracking  target  2  (the  middle  target)  100%  of  the 
time,  and  targets  1  and  3  (the  outer  targets)  less  than  10%  of  the 
time.  (Overall,  ML-PMHT  is  tracking  at  least  one  target  100% 
of  the  time.) 

The  reason  ML-PDA  is  not  performing  as  well  as 
ML-PMHT  in  this  case  is  the  sequential  single-target  im¬ 
plementation  that  is  necessary  for  ML-PDA.  Fig.  17  illustrates 
what  is  happening  in  scenario  2.  This  figure  shows  measure¬ 
ments  from  three  targets  in  close  proximity  to  each  other,  all 
moving  in  the  same  direction.  For  clarity  of  illustration,  the 
measurements  have  no  noise,  and  there  are  no  clutter  mea¬ 
surements  shown.  ML-PDA  is  implemented  in  a  sequential 
mode;  the  first  track  “finds”  the  high-SNR  measurement  in 
each  scan,  and  then  these  measurements  are  excised  from  the 
data.  The  next  track  finds  the  remaining  (relatively)  high-SNR 
measurement,  and  again,  these  measurements  are  excised  from 
the  data.  Since  both  found  tracks  are  using  data  from  all  three 
true  targets,  the  tracks  tend  to  zigzag  back  and  forth  over  the 
middle  target  (and  each  other).  As  a  result,  the  track  scoring 
will  most  likely  result  in  the  middle  target  being  tracked  by 
both  tracks.  Over  many  trials,  this  will  produce  a  high  Pd  value 
for  the  middle  target,  a  duplicate  track  for  the  middle  target, 
and  low  Pd  for  the  outer  targets.  This  is  exactly  what  is  seen 
in  the  Monte  Carlo  results  for  scenario  2:  for  ML-PDA  the 
middle  target  had  a  Pd  of  100%,  while  the  outer  targets  had  a 
Pd  of  less  than  10%.  The  middle  target  also  had,  on  average, 
1.07  duplicate  tracks  associated  with  it. 

In  contrast,  ML-PMHT,  with  its  more  natural  and  appealing 
true  multitarget  formulation,  can  better  find  all  three  targets. 
Simultaneously  optimizing  for  multiple  targets  at  once  prevents 
the  “claiming”  of  the  high-SNR  measurement  by  the  first  track 
to  run  through  the  data;  instead,  the  high-SNR  measurements 
are  more  equally  (and  correctly)  divided  between  the  three 
tracks,  as  is  shown  in  Fig.  18.  Scenario  5  showed  similar  results. 
ML-PDA  was  only  able  to  track  both  targets  approximately 
45%  of  the  time,  while  ML-PMHT  was  able  to  track  both 
targets  approximately  65%  of  the  time.  This  number  is  actually 
slightly  misleading;  the  relatively  low  numbers  for  ML-PDA 


•  High  SNR  measurement  •  Medium  SNR  measurement 
•  Low  SNR  measurement 


Target  Motion  (for  all  three  targets) - ^ 

•- — 

%.  -•ft. 

- -Track  1  - Track  2  .  Track  3 


Fig.  18.  Example  of  multitrack  measurement  assignment  for  ML-PMHT. 

are  because  tracks  using  this  algorithm  switched  targets  on  the 
portion  of  the  trajectory  where  the  targets  were  close,  and  the 
tracks  ended  up  on  the  wrong  targets,  driving  down  Pd.  There 
was  some  switching  for  ML-PMHT  as  well,  but  as  is  seen  by 
the  results,  ML-PDA  was  more  susceptible  to  this  problem. 
This  is  a  result  of  the  sequential  single-target  implementation 
used  for  ML-PDA  that  is  described  above. 

Individual  examples  from  the  Monte  Carlo  runs  again  rein¬ 
force  this  conclusion.  A  scenario  2  result  for  ML-PDA  is  shown 
in  Fig.  7.  Here,  the  first  track  (on  the  middle  target)  “claims”  and 
then  excises  all  the  high-SNR  measurements.  A  second  track 
is  initiated,  and  it  claims  the  remaining  (relatively)  high-SNR 
measurements,  in  the  process  crossing  over  the  first  track  sev¬ 
eral  times.  As  a  result,  it  appears  that  this  track  would  be  as¬ 
sociated  with  the  middle  target  as  well.  There  are  not  enough 
high-SNR  measurements  left  to  form  a  third  track.  In  contrast, 
the  multitarget  ML-PMHT  tracker,  shown  in  Fig.  8,  is  able  to 
track  all  three  targets  with  a  minimum  of  switching.  Examples 
from  scenario  5  show  similar  results.  In  Fig.  13,  ML-PDA  is 
not  able  to  separate  the  targets  on  the  portion  of  the  track  where 
the  two  targets  are  paralleling  each  other  at  very  close  range; 
during  this  portion,  the  tracks  actually  switch  between  targets 
five  times.  In  contrast,  in  Fig.  14,  the  multitarget  ML-PMHT 
tracker  is  able  to  track  both  targets  without  any  switching. 

Finally,  ML-PMHT  outperformed  ML-PDA  on  scenario  4 
for  almost  all  of  the  targets  for  similar  reasons.  (Results  from 
only  targets  7  and  12  are  shown  below;  they  are  fairly  rep¬ 
resentative  of  all  13  targets  in  this  scenario).  Again,  the  se¬ 
quential  single-target  tracking  framework  for  ML-PDA  could 
not  perform  as  well  as  the  true  multitarget  implementation  for 
ML-PMHT.  Tracks  using  the  ML-PDA  implementation  were 
far  more  susceptible  to  being  drawn  off  from  one  target  to  an¬ 
other  by  high-SNR  measurements  (similar  to  the  effect  illus¬ 
trated  in  Fig.  17).  Compare  the  performance  of  ML-PDA  on 
this  scenario  in  Fig.  11  with  the  performance  of  ML-PMHT  in 
Fig.  12.  Several  instances  of  track  switching  are  visible  in  the 
ML-PDA  plot  that  are  not  seen  in  the  ML-PMHT  plot.  While 
this  switching  is  not  between  targets  with  similar  dynamics  (as 
happens  with  scenarios  2  and  5),  it  is  happening  with  targets  that 
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are  moving  very  slowly.  As  a  result,  there  are  no  dynamics — ei¬ 
ther  evolving  position  over  time  or  Doppler — to  differentiate 
measurements  between  two  targets. 

All  the  simulations  up  to  this  point  were  performed  using 
the  ML-PDA  target  measurement  generation  model,  that  is,  at 
most  one  measurement  was  generated  by  the  target  in  any  given 
scan.  In  this  condition,  ML-PMHT  and  ML-PDA  had  identical 
performance  in  the  single-target  cases,  and  ML-PMHT  outper¬ 
formed  ML-PDA  in  the  small-separation  multitarget  cases.  We 
now  round  out  the  comparison  between  the  two  algorithms  by 
considering  what  happens  when  the  ML-PMHT  target  measure¬ 
ment  generation  model  is  used  to  generate  the  data  and  more 
than  one  measurement  is  allowed  to  originate  from  the  target. 

To  do  this,  scenario  1  (with  iT-distributed  clutter)  was  rerun 
for  one,  two,  and  three  expected  target  returns  per  scan.  Re¬ 
sults  for  all  metrics  were  virtually  the  same  between  the  algo¬ 
rithms  (and  similar  to  results  shown  above  for  scenario  1),  with 
the  exception  of  the  number  of  duplicate  tracks.  These  results 
are  shown  in  Table  XII.  Here,  ML-PDA,  as  expected,  suffers 
from  an  increasing  number  of  duplicate  tracks  as  A,  the  ex¬ 
pected  number  of  target  measurements  per  scan,  is  increased. 
In  contrast,  ML-PMHT  basically  has  on  average  no  duplicate 
tracks  for  any  number  of  target  measurements  per  scan.  With 
this  model  of  target-measurement  generation,  ML-PMHT  is  the 
superior  algorithm. 

VI.  Conclusion 

We  developed  a  true  multitarget  implementation  of 
ML-PMHT.  In  the  process,  we  developed  an  expression 
for  the  ML-PMHT  CRLB.  We  showed  that  ML-PMHT  is  a 
statistically  efficient  estimator,  and  that  the  ML-PMHT  CRLB 
can  be  evaluated  in  real  time,  a  valuable  result  in  and  of  itself; 
this  allows  for  a  “complete”  tracker  framework  that  provides 
both  a  state  estimate  as  well  as  a  reliable  estimate  of  the  state 
covariance. 

Then,  we  tested  these  developments  by  comparing  ML-PDA 
and  ML-PMHT  with  Monte  Carlo  testing.  This  testing  first 
showed  that  ML-PDA  and  ML-PMHT  are  effective  VLO 
trackers,  working  down  to  an  expected  target  SNR  of  4-5 
dB  (postsignal  processing).  After  this,  the  ML-PDA  and 
ML-PMHT  tracking  algorithms  were  applied  to  five  different 
benchmark  scenarios  with  Monte  Carlo  trials  using  a  target 
measurement  generation  model  of  zero  or  one  measurements 
originating  from  the  target  in  a  scan.  For  scenarios  with  a  single 
target  or  multiple  targets  with  measurements  that  could  easily 
be  differentiated  by  dynamics,  the  performances  of  ML-PDA 
and  ML-PMHT  were  identical:  ML-PMHT  did  not  suffer  from 
the  fact  that  its  measurement  assignment  model  did  not  match 
the  actual  target  measurement  generation  model.  In  cases  with 
closely  spaced  targets  with  measurements  that  could  not  be 
differentiated  easily  by  dynamics,  ML-PMHT  outperformed 
ML-PDA  because  the  former  had  a  true  multitarget  LLR 


formulation,  while  the  latter  had  to  handle  multiple  targets  in 
a  sequential  single-target  mode.  Finally,  when  the  target  mea¬ 
surement  generation  model  was  switched  to  that  of  ML-PMHT, 
with  multiple  measurements  per  scan  being  generated  by  the 
target,  ML-PMHT  outperformed  ML-PDA  in  terms  of  the 
number  of  duplicate  tracks  generated.  Overall,  the  performance 
of  ML-PMHT  makes  it  the  preferred  algorithm. 


Appendix 

Calculation  of  the  ML-PMHT  FIM 

Here,  we  fully  derive  the  FIM  for  a  window  of  ML-PMHT 
data.  This  is  a  specific  implementation  of  the  general  work  given 
in  [22].  The  derivation  is  also  similar  to  that  done  for  ML-PDA 
in  [8].  However,  the  ML-PMHT  calculation  ends  up  being  at 
most  a  fourfold  integral  (as  opposed  to  approximately  an  ra-fold 
integral  for  ML-PDA,  where  m  is  the  number  of  measurements 
in  a  scan),  so  this  can  be  done  in  real  time  instead  of  having 
to  calculate  it  offline  with  Monte  Carlo  integration,  as  must  be 
done  with  ML-PDA. 

Consider  a  window  of  ML-PMHT  data  with  Nw  scans.  Since 
all  measurements  from  scan  to  scan  are  assumed  to  be  indepen¬ 
dent,  the  FIM  J  of  the  total  window  can  just  be  written  as  the 
sum  of  the  FIM’s  J;  of  the  individual  scans 

Nw 

j  =  X>-  <29> 

Now,  the  FIM  of  a  scan  is  calculated  as 
J i  =  E{[Vx  In p[Z(i)  |  x]][Vx  In p[Z(i)  |  x]]T}|x=Xtrue.  (30) 


The  ML-PMHT  likelihood  for  a  single  scan  is 

m  i 

v\Z{i)  |x]  =  n  {yPoKO  +  7TlP[Zj(j)  |x]p[(o,-)}  . 

3  =  1 

Taking  the  gradient  of  the  logarithm  of  this  gives 

m7; 

Vxlnp[Z(i)\x]  =  Y/V^ 

.7  =  1 

x  VEM _ J 

KoPofaj)  +  TTl Pl(flj)  p-l( zJ-^)rR71(zJ-^)  ' 

V 


(31) 


(32) 


Here,  (j)  is  the  state-to-measurement  conversion  and  is  the 
Jacobian  of  0.  (Note  that  we  have  dropped  the  scan  index  i 
on  the  values  z,  </>,  R,  and  a  for  ease  of  expression.)  Inserting 
(32)  into  (30)  and  recognizing  that  cross  terms  in  the  product 
of  the  sums  will  go  to  zero  results  in  (33),  shown  at  the  bottom 
of  the  page,  where  dZ  =  n[=i  dzi  and  da.  =  YIi=i  dai .  We 
are  denoting  the  denominator  in  (32)  as  p(zj  |  x).  Note  that  all 
p(zi  |  x)  terms  for  l  /  j  are  separable  and  simply  integrate  to 


KiPlUj)] 


J*  =  X 


1(«J-^)R-1  (Zj  _  ^)(z>  _  ^.)tR-1 

p(  zi  lx)2 


-B^Y[p(zi  |  x)dZda 


(33) 


316 


17 


IEEE  JOURNAL  OF  OCEANIC  ENGINEERING,  VOL.  39,  NO.  2,  APRIL  2014 


rrii  r 

=  E/D 

.7  =  1  JV 


^oPoUj)  +  ^l/U/'i),,  ‘.(v,,  <f>) '  k  • 

v  fW^I\ 


-'D^dzjdaj 


(34) 


one,  and  the  term  for  l  —  j  is  canceled  by  one  of  the  terms  in  the 
denominator.  This  happens  because  the  likelihood  function  for 
ML-PMHT  (for  a  scan)  is  a  product,  and  this  reduces  the  overall 
FIM  integral  from  being  an  ra-fold  integral  to  only  a  threefold 
or  fourfold  integral.  In  contrast,  as  is  shown  in  [8],  the  ML-PDA 
likelihood  function  for  a  scan  of  data  consists  of  a  summation, 
not  a  product,  over  all  the  measurements  in  the  scan.  As  a  result, 
taking  the  gradient  of  the  logarithm  of  this  results  in  summation 
over  all  the  measurements  in  the  denominator  of  the  ML-PDA 
equivalent  of  (33).  Thus,  the  cancellation  that  occurs  above  does 
not  happen  for  ML-PDA;  its  FIM  integral  remains  on  the  order 
of  m-fold. 

At  this  point,  we  can  reduce  (33)  to  (34),  shown  at  the  top 
of  the  page.  The  Jacobians  are  independent  of  the  variables  of 
integration  and  the  summation  indices,  and  thus  can  be  moved 
outside  the  integral  and  summation.  Furthermore,  we  can  make 
the  following  substitution  of: 


0  =  G;)(z(  t)  (35) 

where  Gj  is  the  Cholesky  decomposition  of  R J1  such  that 
Gj  Gj  =  R71.  With  this,  it  is  possible  to  write  the  final  ex¬ 
pression  for  J  i 

m  i 

j*  =  d£EgJ 

3  =  1 
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(36) 


This  is  a  relatively  simple  threefold  or  fourfold  integration  (de¬ 
pending  on  the  dimension  of  z j)  in  this  application  and  can 
be  done  in  real  time.  The  resultant  J;  for  each  scan  is  simply 
summed  up  then  in  accordance  with  (29)  to  get  the  FIM  for  a 
window  of  data. 
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